Swap_Mult Analysis

This is an R Markdown Notebook. When you execute code within the notebook, the results appear beneath the code.

# removes all variables but NOT functions
rm(list = setdiff(ls(), lsf.str()))

knitr::opts_chunk$set(echo = FALSE)
library(magrittr)
library(ggplot2)
library(dplyr)

Attaching package: ‘dplyr’

The following objects are masked from ‘package:stats’:

    filter, lag

The following objects are masked from ‘package:base’:

    intersect, setdiff, setequal, union

library(ggplot2)
library(lme4)

Loading required package: Matrix

library(tidyr)

Attaching package: ‘tidyr’

The following object is masked from ‘package:Matrix’:

    expand

The following object is masked from ‘package:magrittr’:

    extract

library(merTools) 

Loading required package: arm
Loading required package: MASS

Attaching package: ‘MASS’

The following object is masked from ‘package:dplyr’:

    select


arm (Version 1.9-3, built: 2016-11-21)

Working directory is /Users/djw/Dropbox/PROGRAMMING/_NEURO/2017_MADE/Analysis

library(lsr)
library(reshape2)

Attaching package: ‘reshape2’

The following object is masked from ‘package:tidyr’:

    smiths

library(ggpubr)
library(plyr)

--------------------------------------------------------------------------------------------------------------------
You have loaded plyr after dplyr - this is likely to cause problems.
If you need functions from both plyr and dplyr, please load plyr first, then dplyr:
library(plyr); library(dplyr)
--------------------------------------------------------------------------------------------------------------------

Attaching package: ‘plyr’

The following objects are masked from ‘package:dplyr’:

    arrange, count, desc, failwith, id, mutate, rename, summarise, summarize

library(Hmisc) # cut2

Loading required package: lattice
Loading required package: survival
Loading required package: Formula

Attaching package: ‘Hmisc’

The following objects are masked from ‘package:plyr’:

    is.discrete, summarize

The following objects are masked from ‘package:dplyr’:

    combine, src, summarize

The following objects are masked from ‘package:base’:

    format.pval, round.POSIXt, trunc.POSIXt, units

library(RColorBrewer)

specify_decimal <- function(x, k) trimws(format(round(x, k), nsmall=k))

Image Swap Version

Analysis of the MULTIPLIER trials first

still need to work out how to select subject if the NM trials based on their performance in the multiplier trials

Dataframe header

# set WD
setwd("~/Dropbox/PROGRAMMING/_NEURO/2017_MADE/Analysis")
load("Data/S_M.Rdata")
load("Data/S_M_raw.Rdata")
load("Data/NS_M.Rdata")
load("Data/S_NM.Rdata")
load("Data/S_NM_raw.Rdata")
head(S_M_raw)

mean(S_NM$rt)
mean(S_NM$accuracy)
mean(S_NM$swapAmount)

mean(S_M$rt)
mean(S_M$accuracy)
sd(S_M$accuracy)
mean(S_M$swapAmount)


subject_means <- group_by(S_M_raw, subject) %>%
  dplyr::summarize(swaps = mean(swapAmount, na.rm =T))
plot(subject_means)

subject_means <- group_by(S_M_raw, subject) %>%
  dplyr::summarize(rt = mean(rt, na.rm = T), finalEarnings = mean(finalEarnings, na.rm = T))

Table of subject count, mean, median, SD, range, skew, kurtosis

Experiment 1/2

Accuracy, rt

Uncleaned/Cleaned

load("Data/made_descripive_stats.Rdata")

Add all DFs (cleaned), mark by condition

Plot boxplots of rt and accuracy

load("Data/NS_NM.Rdata")
load("Data/NS_M.Rdata")
load("Data/S_NM.Rdata")
load("Data/S_M.Rdata")
NS_NM$type = "NoSwap_NoMult"
NS_M$type = "NoSwap_Mult"
S_NM$type = "Swap_NoMult"
S_M$type = "Swap_Mult"
# Join all DFs
common <- intersect(names(NS_NM), names(NS_M))
df = rbind(NS_NM[,common], NS_M[,common])
common <- intersect(names(S_NM), names(S_M))
df2 = rbind(S_NM[,common], S_M[,common])
common <- intersect(names(df), names(df2))
df_full = rbind(df[,common], df2[,common])
# Boxplot RT
boxplot(rt ~ factor(type),
        varwidth = TRUE, xlab = "Trial Type",
        main = "RT by Trial Type", ylab = "RT in s", data = df_full)

# Barplot of Accuracy
#First get means for each trial condition (type) by Subject
d <- df_full
subject_means <- group_by(d, type) %>%
  dplyr::summarize(accuracy = mean(correct, na.rm = T))
subject_means$accuracy = round(subject_means$accuracy, 3)
#PLOT
lower=c(0.787990, 0.787694, 0.798404, 0.759005) 
upper=c(0.877638, 0.867699, 0.880029, 0.893729)
barplot <- ggplot(subject_means, aes(x = type, y = accuracy)) +
  stat_summary(
    geom = "bar",
    fun.y = "mean",
    col = "black",
    fill = "gray70"
  ) +
  #geom_point(position = position_jitter(h = 0, w = 0.2)) +
  scale_y_continuous(limits = c(0, max(d$correct, na.rm = T)),
                     expand = c(0, 0)) +
  geom_text(aes(label=accuracy), position=position_dodge(width=0.9), vjust=-0.25, hjust=-0.65) +
  geom_errorbar(data=subject_means, mapping=aes(x=type, ymin=upper, ymax=lower), width = 0.3)
barplot + ggtitle("Accuracy by Trial Type")

T-Test NS_M vs S_M accuracy

t.test(NS_M$correct, S_M$correct)

t.test(NS_NM$rt[NS_NM$Trial<51], NS_NM$rt[NS_NM$Trial>50])

subject_means <- group_by(S_M_raw, subject) %>%
  dplyr::summarize(swapAvg = mean(swapAvg), finalEarnings = mean(finalEarnings, na.rm = T))
x = subject_means$swapAvg[subject_means$swapAvg>1.6]
length(x)
length(subject_means$subject)

T-test NS vs S RT

t.test(subject_means$rt[subject_means$study == "Standard Mult."], subject_means$rt[subject_means$study=="Swap Mult."])

    Welch Two Sample t-test

data:  subject_means$rt[subject_means$study == "Standard Mult."] and subject_means$rt[subject_means$study == "Swap Mult."]
t = -2.8654, df = 35.986, p-value = 0.006914
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -1.1442602 -0.1957772
sample estimates:
mean of x mean of y 
 2.361594  3.031612 

There were 50 or more warnings (use warnings() to see the first 50)

Summed Val vs. Accuracy

datac <- summarySEwithin(df, measurevar="correct", withinvars=c("multDif","difficulty"), idvar="subject")

Automatically converting the following non-factors to factors: multDif, difficulty

ggplot(datac, aes(x=difficulty, y=correct, fill=multDif)) +
    geom_bar(position=position_dodge(.9), colour="black", stat="identity") +
    geom_errorbar(position=position_dodge(.9), width=.25, aes(ymin=correct-ci, ymax=correct+ci)) +
    coord_cartesian(ylim=c(0.0,1)) +
    labs(y = "Accuracy", x = "Difficulty (net value)") +
    scale_x_discrete(labels=c("1" = "Easy (>=0.5)", "2" = "Difficult (<0.5)")) +    
    scale_y_continuous(breaks=seq(0,1,0.1)) +
    theme_bw() +
    theme(axis.title.x=element_text(size=18),
        axis.title.y = element_text(size = 18)) +
    scale_fill_brewer(palette="Pastel2") +
    guides(fill=guide_legend(title="Difference\nbetween\nMultipliers"))
setwd("/Users/djw/Dropbox/PHD/PRESENTATIONS/2017_SNE/Plots")

The working directory was changed to /Users/djw/Dropbox/PHD/PRESENTATIONS/2017_SNE/Plots inside a notebook chunk. The working directory will be reset when the chunk is finished running. Use the knitr root.dir option in the setup chunk to change the the working directory for notebook chunks.

ggsave("MultDif.pdf", width = 16, height = 12, units = "cm")

Summed val vs. RT

df <- S_M

There were 50 or more warnings (use warnings() to see the first 50)

#RT vs. Summed Value
ggplot() +
  geom_smooth(aes(x=summedVal, y=rt, group = factor(multNum), colour = factor(multNum)), df) +
  #geom_smooth(aes(x=summedVal, y=rt, colour = "flip"), subset(df, flip==1)) +
  coord_cartesian(xlim = c(-3, 3))  +
  #ggtitle("RT vs. Summed Val")
  #geom_point(shape=1) +    # Use hollow circles
  geom_smooth() +  # Add a loess smoothed fit curve with confidence region
  theme_minimal()+
  guides(colour=guide_legend("Multiplier \nCondition")) +
  scale_x_continuous(name="Net Value ($)", seq(-3,3,0.5), limits = c(-3,3))+
  scale_y_continuous(name = "Reaction Time (s)")

setwd("/Users/djw/Dropbox/PHD/PRESENTATIONS/2017_SNE/Plots/")

The working directory was changed to /Users/djw/Dropbox/PHD/PRESENTATIONS/2017_SNE/Plots inside a notebook chunk. The working directory will be reset when the chunk is finished running. Use the knitr root.dir option in the setup chunk to change the the working directory for notebook chunks.

ggsave("2ndFix.pdf", width = 18, height = 12, units = "cm")

Second Fixation vs Summed Value

df <- S_M
#RT vs. Summed Value
ggplot() +
  geom_smooth(aes(x=summedVal, y=`2_fixation`, group = factor(secondMult), colour = factor(secondMult)), df) +
  #geom_smooth(aes(x=summedVal, y=logRT, colour = "flip"), subset(total_M_clean3, flip==1)) +
  coord_cartesian(xlim = c(-3, 3))  +
  #ggtitle("Second Fixation vs Summed Value")
  #geom_point(shape=1) +    # Use hollow circles
  geom_smooth()  +# Add a loess smoothed fit curve with confidence region
  theme_minimal()+
  guides(colour=guide_legend("Multiplier \nCondition")) +
  scale_x_continuous(name="Trial Net Value ($)", seq(-3,3,0.5), limits = c(-3,3))+
  scale_y_continuous(name = "Fixation Duration (s)")

#Test for SIG
summary(lm(`2_fixation`~summedVal + factor(secondMult), df))

Call:
lm(formula = `2_fixation` ~ summedVal + factor(secondMult), data = df)

Residuals:
    Min      1Q  Median      3Q     Max 
-1.0029 -0.3556 -0.1513  0.1609  6.7574 

Coefficients:
                     Estimate Std. Error t value Pr(>|t|)    
(Intercept)          1.034731   0.009812 105.452  < 2e-16 ***
summedVal           -0.012915   0.006777  -1.906   0.0567 .  
factor(secondMult)2  0.115010   0.019281   5.965 2.57e-09 ***
factor(secondMult)3  0.184020   0.019372   9.499  < 2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 0.6114 on 6576 degrees of freedom
  (17 observations deleted due to missingness)
Multiple R-squared:  0.01601,   Adjusted R-squared:  0.01556 
F-statistic: 35.67 on 3 and 6576 DF,  p-value: < 2.2e-16

QUESTIONNAIRES

Boxplot of unfiltered data

You can see which people were not trying:

df <- S_M_raw

# Filter out extremely long rt times
df <- df[!(df$rt>10),]

boxplot(rt ~ factor(subject),
        varwidth = TRUE, xlab = "subject",
        main = "Boxplot of RT conditional on\
        subject", ylab = "RT", data = df)

Plotting mean RT vs. Earnings (unfiltered data):

KMeans to Divide into groups based on swaps/$$ (unfiltered data):

Tends to include some people who earned over $50 but never looked at the second image…

Clean data by removing people with less than 1.5 swap average and less than $0 earnings

#remove earnings below 0
total_M_clean2 <- S_M_raw[!(S_M_raw$finalEarnings<0),]
#remove flip avgs less than 1.5
total_M_clean2 <- total_M_clean2[!(total_M_clean2$flipAvg<1.5),]
length(unique(total_M_clean2$subject))

#MORE AGRESSIVE: REMOVE BELOW 75%
total_M_clean3 <- S_M_raw[!(S_M_raw$accuracy<0.75),]
length(unique(total_M_clean3$subject))

Taking out those that swapped images less than 2 times and earned less than $0 leaves 24 subjects Taking out those that had an accuracy of below 75% leaves 21 subjects

Boxplot of filtered data

A lot of variation in avg. RT and the variance of RTs. 4,7 and 23 seem to have particularly low variance.

boxplot(rt ~ factor(subject),
        varwidth = TRUE, xlab = "subject",
        main = "Boxplot of RT conditional on\
        subject", ylab = "RT", data = S_M)

Check RT distribution

Skewed right, as is typical of RT.

Histogram of

There were 41 warnings (use warnings() to see them)

is.na() applied to non-(list or vector) of type 'NULL'is.na() applied to non-(list or vector) of type 'NULL'

RTs By Subject

g = ggplot(S_M, aes(x = rt)) + geom_histogram() 

g + facet_wrap(~ subject, ncol=2)

Hist of Middle Fixations (All Subjects)

df <- S_M

#create version of dataset with only middle fixations (we use this later for drawing random middle fixations)

# make sure every trial ends with an NA val
df$`14_fixation` <- NA     

# find col Number of first fix
# Check that we are starting at column fix#1
firstFix = which(names(df) == "1_fixation")

# find the first NA val in each row
naVal <- vector(mode="numeric", length=0)
for(x in 1:length(df$Trial)){
  naVal[x] <- min(which(is.na(df[x,])))
}

# if naVal >39 then there were at least 3 fixations
midFix <- vector(mode = "numeric", length = 0)
for(x in 1:length(df$Trial)){
  if(naVal[x]>(firstFix+2)){
    for(y in (firstFix+1):(naVal[x]-2)){
      val = df[x,y]
      midFix <- c(midFix, val)
    }
  }
}

mean(midFix)
median(midFix)

S_M$finalFix
hist(midFix,
   col=rgb(1,0,0,0.5), breaks=seq(0,9,0.1), ylim=c(0,900), xlab="RT", main = "Middle Fixations")
#abline(v=median(midFix), col="blue")

legend("topright",
       c(as.expression(bquote(MidFix_Median == .(median(midFix))))),
       col = c("blue"),
       lwd = c(2, 2, 1))

Hist of Middle Fixations by Subject

load("Data/S_M_K.Rdata")
df <- S_M_K

# REMOVE First and Last Fixations
df <- df[(df$fixNum!=1 & df$revFixNum!=1),]

# SUBJECT MEANS AND SD AND FIX #
subject_means <- group_by(df, subject) %>%
  dplyr::summarize(median = median(rt, na.rm = T), mean = mean(rt, na.rm = T), sd = sd(rt, na.rm = T), count = length(rt))
subject_means

# Plot Hist
g = ggplot(df, aes(x = rt)) + 
  geom_histogram(data = transform(df, subject = NULL), fill = "blue", alpha = 0.4) +
  geom_histogram() 

g + facet_wrap(~ subject, ncol=2)

# PLot Density
# "adjust" controls the bandwidth
allFix = transform(df, subject = NULL)

g = ggplot(df, aes(x=rt)) + 
  geom_density(data = allFix, aes(alpha = 0.5, fill = "group")) +
  geom_density(aes(alpha = 0.5, fill = "subject")) +
  #geom_vline(aes(xintercept=mean(rt)), color="black", size=1) +
  scale_fill_manual(name = "Density Plot", 
                    values = c(group = "blue", subject = "red"))

appender <- function(string, prefix = "Subject: ", suffix = "  Mean: ", mean = specify_decimal(subject_means$mean[subject_means$subject == (as.numeric(string))],2)) paste0(prefix, string, suffix, mean)
              
g + facet_wrap(~ subject, ncol=2, labeller = as_labeller(appender)) +
  theme(strip.text = element_text(size = 22),
        axis.text.x = element_text(size = 18),
        axis.title.x = element_text(size = 25),
        axis.text.y = element_text(size = 18),
        axis.title.y = element_text(size = 25),
        legend.text = element_text(size = 18),
        legend.title = element_text(size = 21)) +
  ggtitle("Subject vs Group Middle Fixations (in s)")

Log RT histogram

Plot RT vs Summed Value

Plot LOG RT vs Summed Value

Plot Accuracy vs Difficulty

Hist of Accuracy vs Difficulty

SIGMOID of ACCEPT VS VALUE

d <- S_M

plot(d$summedVal, d$choice,
     main = "Choice vs. Summed Val",
     xlab="Summed Val", ylab="P (Accept)",
     xlim=c(-3, 3))

model <- glm(choice ~ summedVal, data=d, family=binomial(link = logit))
summary(model)

xv <- seq(min(d$summedVal), max(d$summedVal), 0.01)
yv <- predict(model,list(summedVal=xv), type="response")

abline(0.5,0, lty=2)
lines(xv,yv)

#Find the inflection point where there is a 50/50 probability of subject accepting.
p <- 0.5
x <- (log(p/(1-p)) - coef(model)[1]) / coef(model)[2]
x

Final Fix vs Choice

TWO DIFFICULTY CONDITIONS (Hard/Easy)

Plot Hard vs. Easy for Subjects

Plot LOG RT vs Absolute Summed Value

Plot Image Swaps vs Summed Value

Plot Image Swaps vs Absolute Value

Mean Number of Fixations

Summed Value vs. Accept/Reject Decision

Logistic regression curve

RT Correct vs Incorrect Responses

#GOOD BASE FOR PLOTTING EXAMPLE

#BAR PLOT
d <- S_M
d$correct[d$correct==1] = "Correct"
d$correct[d$correct==0] = "Incorrect"

subject_means <- group_by(d, subject, correct) %>%
  dplyr::summarize(rt = mean(rt, na.rm = T))
subject_means

#PLOT
barplot <- ggplot(subject_means, aes(x = correct, y = rt, fill=correct)) +
  stat_summary(
    geom = "bar",
    fun.y = "mean",
    col = "black",
    stat = "identity"
  ) +
  geom_point(position = position_jitter(h = 0, w = 0.2)) +
  scale_y_continuous(limits = c(0, max(d$rt, na.rm = T)),
                     expand = c(0, 0))+
  labs(y = "RT (seconds)", x = "Response")+
  theme_minimal()+
  theme(legend.position="none") +
  ggtitle("RT vs Correct Choice")

barplot 
t.test(d$rt[d$correct=="Correct"], d$rt[d$correct=="Incorrect"])

BarPlot for Choice vs RT

#BAR PLOT
d<- S_M
d$choice[d$choice==1] = "Accept"
d$choice[d$choice==0] = "Reject"

subject_means <- group_by(d, subject, choice) %>%
  dplyr::summarize(rt = mean(rt, na.rm = T))
subject_means

#PLOT
barplot <- ggplot(subject_means, aes(x = choice, y = rt, fill=choice)) +
  stat_summary(
    geom = "bar",
    fun.y = "mean",
    col = "black",
    stat = "identity"
  ) +
  geom_point(position = position_jitter(h = 0, w = 0.2)) +
  scale_y_continuous(limits = c(0, max(d$rt, na.rm = T)),
                     expand = c(0, 0))+
  labs(y = "RT (seconds)", x = "Response") +
  theme_minimal() +
  theme(legend.position="none") +
  ggtitle("RT vs Choice")
  #stat_compare_means(label.y = 8.0) +
  #stat_compare_means(ref.group = "Accept", label = "p.signif", label.y = c(7.0))

barplot 

t.test(d$rt[d$choice=="Accept"], d$rt[d$choice=="Reject"])

Deciding Factor: House/Face

#BAR PLOT
d<- S_M
# Reduce DF to Decider Trials
x<- (d$faceTotal>0 & d$houseTotal<0) | (d$faceTotal<0 & d$houseTotal>0)
d = d[x==TRUE,]
# Create Decider Column
d$decider = 0
d$decider[abs(d$faceTotal)>abs(d$houseTotal)] = "Face"
d$decider[abs(d$houseTotal)>abs(d$faceTotal)] = "House"
#remove values of zero
d = d[d$decider!=0,]

subject_means <- group_by(d, subject, decider) %>%
  dplyr::summarize(rt = mean(rt, na.rm = T))
subject_means

#PLOT
barplot <- ggplot(subject_means, aes(x = decider, y = rt, fill=decider)) +
  stat_summary(
    geom = "bar",
    fun.y = "mean",
    col = "black",
    stat = "identity"
  ) +
  geom_point(position = position_jitter(h = 0, w = 0.2)) +
  scale_y_continuous(limits = c(0, max(d$rt, na.rm = T)),
                     expand = c(0, 0))+
  labs(y = "RT (seconds)", x = "Response") +
  theme_minimal() +
  theme(legend.position="none") +
  ggtitle("RT vs Decider Stimulus")
  #stat_compare_means(label.y = 8.0) +
  #stat_compare_means(ref.group = "Accept", label = "p.signif", label.y = c(7.0))

barplot 

t.test(d$rt[d$decider=="Face"], d$rt[d$decider=="House"])

Accuracy vs Neg/Pos Summed Value

#BAR PLOT
d<- S_M
# Reduce DF to Decider Trials
d$posNeg <- 0
d$posNeg[d$summedVal>0] = "Positive"
d$posNeg[d$summedVal<0] = "Negative"
d <- d[d$summedVal!=0,] # remove values of 0

subject_means <- group_by(d, subject, posNeg) %>%
  dplyr::summarize(accuracy = mean(correct, na.rm = T))
subject_means
subject_means$meanAcc = mean(subject_means$accuracy)

#PLOT
barplot <- ggplot(subject_means, aes(x = posNeg, y = accuracy, fill=posNeg, label=accuracy)) +
  stat_summary(
    geom = "bar",
    fun.y = "mean",
    col = "black",
    stat = "identity"
  ) +
  geom_point(position = position_jitter(h = 0, w = 0.2)) +
  scale_y_continuous(limits = c(0, max(d$accuracy+0.25, na.rm = T)),
                     expand = c(0, 0))+
  labs(y = "Accuracy", x = "Summed Value") +
  theme_minimal() +
  theme(legend.position="none") +
  ggtitle("Accuracy vs Summed Value Sign")
  #stat_compare_means(label.y = 8.0) +
  #stat_compare_means(ref.group = "Accept", label = "p.signif", label.y = c(7.0))

barplot 

t.test(d$correct[d$posNeg=="Positive"], d$correct[d$posNeg=="Negative"])

RT vs Positive/Negative Summed Value

#BAR PLOT
d<- S_M
# Reduce DF to Decider Trials
d$posNeg <- 0
d$posNeg[d$summedVal>0] = "Positive"
d$posNeg[d$summedVal<0] = "Negative"
d <- d[d$summedVal!=0,] # remove values of 0

subject_means <- group_by(d, subject, posNeg) %>%
  dplyr::summarize(rt = mean(rt, na.rm = T))
subject_means

#PLOT
barplot <- ggplot(subject_means, aes(x = posNeg, y = rt, fill=posNeg)) +
  stat_summary(
    geom = "bar",
    fun.y = "mean",
    col = "black",
    stat = "identity"
  ) +
  geom_point(position = position_jitter(h = 0, w = 0.2)) +
  scale_y_continuous(limits = c(0, max(d$rt-2.5, na.rm = T)),
                     expand = c(0, 0))+
  labs(y = "RT (seconds)", x = "Summed Value") +
  theme_minimal() +
  theme(legend.position="none") +
  ggtitle("RT vs Summed Value Sign")
  #stat_compare_means(label.y = 8.0) +
  #stat_compare_means(ref.group = "Accept", label = "p.signif", label.y = c(7.0))

barplot 

t.test(d$rt[d$posNeg=="Positive"], d$rt[d$posNeg=="Negative"])

First/Middle/Last Fixation mean duration (boxPlot)

load("/Users/djw/Dropbox/PROGRAMMING/*NEURO/aDDM_Krajbich/S_M_K.Rdata")
d <- S_M_K

d$fixType <- 2
d$fixType[d$fixNum==1] <- 1
d$fixType[d$revFixNum==1] <-3

median(d$fixDur[d$fixType==1])
median(d$fixDur[d$fixType==2])
median(d$fixDur[d$fixType==3])


my_comparisons <- list( c("1", "2"), c("1", "3"), c("2", "3") )

p0 = ggboxplot(d, x = "fixType", y = "fixDur", color = "fixType", palette = "jco" )+ #outlier.shape=NA
  
  labs(y = "Fixation Time (seconds)", x = "Fixation Type") +
  theme(legend.position="none") +  
  scale_x_discrete(labels=c("1" = "First", "2" = "Middle",
                              "3" = "Last")) +
  #scale_y_continuous(limits = quantile(d$fixDur, c(0.1, 0.9))) +
  stat_compare_means(comparisons = my_comparisons)+ # Add pairwise comparisons p-value
  stat_compare_means(label.y = 13) +     # Add global p-value
  ggtitle("Fixation Durations")

p0

## How many trials with more than 2 fixatinos
d = S_M
length(d$Trial[d$swapAmount>2])/length(d$Trial)

mean(d$firstVal[d$swapAmount<2])

Summed Val vs Fixations

df <- S_M

#RT vs. Summed Value
ggplot() +
  geom_smooth(aes(x=summedVal, y=swapAmount), df) +
  #geom_smooth(aes(x=summedVal, y=logRT, colour = "flip"), subset(total_M_clean3, flip==1)) +
  coord_cartesian(xlim = c(-3, 3))  +
  ggtitle("Summed Value vs Fixations")
  #geom_point(shape=1) +    # Use hollow circles
  geom_smooth()  # Add a loess smoothed fit curve with confidence region

#Test for SIG
summary(lm(`2_fixation`~summedVal + factor(secondMult), df))

# create a dummy data frame with outliers
df = data.frame(y = c(-100, rnorm(100), 100))

# create boxplot that includes outliers
p0 = ggplot(df, aes(y = y)) + geom_boxplot(aes(x = factor(1)))

p0
# compute lower and upper whiskers
ylim1 = boxplot.stats(df$y)$stats[c(1, 5)]

# scale y limits based on ylim1
p1 = p0 + coord_cartesian(ylim = ylim1*1.05)
p1

Plot RT with/without Flip

T Test on RT difference between flip and non flip trials RT

library(tidyr)
subject_means_wide <-
  spread(subject_means,
         key = flip,
         value = rt,
         sep = "_")
subject_means_wide

#T-TEST for flip vs. non-flip trials
t.test(subject_means_wide$flip_1, subject_means_wide$flip_2, paired = TRUE)
mean(subject_means_wide$flip_2)
sd(subject_means_wide$flip_2)

Plot performance with/without Flip

There were 50 or more warnings (use warnings() to see the first 50)

Based on the plot it looks like flip trials on average do WORSE and have higher VARIANCE between subjects (in addition to taking longer)

T Test on RT difference between flip and non flip trials % Correct

subject_means_wide <-
  spread(subject_means,
         key = flip,
         value = corPct,
         sep = "_")
subject_means_wide

#T-TEST for flip vs. non-flip trials
t.test(subject_means_wide$flip_1, subject_means_wide$flip_2, paired = TRUE)
sd(subject_means_wide$flip_2)

It turns out that the difference in performance is not significant.

Plot RT vs negative/positive summed Val

Based on the plot it looks like RT is longer for negative summed values.

T Test on RT difference between pos/neg summed values

subject_means_wide <-
  spread(subject_means,
         key = posNegSum,
         value = rt,
         sep = "_")

#T-TEST for flip vs. non-flip trials
t.test(subject_means_wide$posNegSum_0, subject_means_wide$posNegSum_1, paired = TRUE)

This suggests that there is a significantly longer time spent on choices with a negative summed value.

Plot % Correct vs negative/positive summed Val

Based on the plot it looks like people do better when the summed val is positive.

T Test on RT difference between flip and non flip trials % Correct

subject_means_wide <-
  spread(subject_means,
         key = posNegSum,
         value = corPct,
         sep = "_")
subject_means_wide

#T-TEST for flip vs. non-flip trials
t.test(subject_means_wide$posNegSum_0, subject_means_wide$posNegSum_1, paired = TRUE)

The difference is significant. So people take less time but perform better when the summed value is positive.

Multipliers and Absolute Net Value

d <- S_M

# EFFECTS on Abs Val due to MULTS
subject_means <- group_by(d, subject, multNum) %>%
  dplyr::summarize(absNet = mean(absSummedVal, na.rm = T), rt = mean(rt, na.rm = T))
subject_means

# Mean by Mult
mean(subject_means$absNet[subject_means$multNum == 0])
mean(subject_means$absNet[subject_means$multNum == 1])
mean(subject_means$absNet[subject_means$multNum == 2])

# SD by Mult
sd(subject_means$absNet[subject_means$multNum == 0])
sd(subject_means$absNet[subject_means$multNum == 1])
sd(subject_means$absNet[subject_means$multNum == 2])

# 0 1
t.test(subject_means$absNet[subject_means$multNum == 0],
                            subject_means$absNet[subject_means$multNum == 1], paired = TRUE)
# 1 2
t.test(subject_means$absNet[subject_means$multNum == 1],
                            subject_means$absNet[subject_means$multNum == 2], paired = TRUE)
# 0 2
t.test(subject_means$absNet[subject_means$multNum == 0],
       subject_means$absNet[subject_means$multNum == 2], paired = TRUE)

RT distributions for different multNums

d = S_M

med.fac = ddply(d, .(multNum), function(.d)
data.frame(x=median(.d$rt)))

# HISTOGRAM VERSION
p <- ggplot(data = d, aes(x = rt, fill=multNum)) + 
  geom_histogram() +
  labs(title="RT Distribution vs. Number of Multipliers", x = "RT (seconds)", y ="Count") +
  geom_vline(data=med.fac, aes(xintercept=x)) +
  theme_minimal() +
  theme(legend.position="none")
p + facet_wrap(~multNum)

# BOXPLOT VERSION

# create medians to insert as text
x <- d
x$multNum = d$multNum+1
p_meds <- ddply(x, .(multNum), summarise, med = median(rt))
p_meds$med = round(p_meds$med, digits = 2)  # round to two decimal values

# List of conditions to compare
my_comparisons <- list( c("0", "1"), c("0", "2"), c("1", "2") )

p0 = ggboxplot(d, x = "multNum", y = "rt", color = "multNum", palette = "jco" )+ #outlier.shape=NA
  
  labs(y = "Total RT (seconds)", x = "Number of Multipliers") +
  theme(legend.position="none") +  
  scale_x_discrete(labels=c("0" = "None", "1" = "One",
                              "2" = "Two")) +
  #scale_y_continuous(limits = quantile(d$fixDur, c(0.1, 0.9))) +
  stat_compare_means(comparisons = my_comparisons)+ # Add pairwise comparisons p-value
  stat_compare_means(label.y = 13) +     # Add global p-value
  geom_text(data = p_meds, aes(x = multNum, y = med, label = med), 
              size = 3, vjust = -1.5) +
  ggtitle("RT Distribution vs. Number of Multipliers")

p0

RT by Mult for “Difficult” trials

# create medians to insert as text
d<- S_M
d <- d[d$absSummedVal<0.5, ] # limit to absolute summed values below 0.50
x <- d
x$multNum = d$multNum+1
p_meds <- ddply(x, .(multNum), summarise, med = median(rt))
p_meds$med = round(p_meds$med, digits = 2)  # round to two decimal values

# List of conditions to compare
my_comparisons <- list( c("0", "1"), c("0", "2"), c("1", "2") )

p0 = ggboxplot(d, x = "multNum", y = "rt", color = "multNum", palette = "jco" )+ #outlier.shape=NA
  
  labs(y = "Total RT (seconds)", x = "Number of Multipliers") +
  theme(legend.position="none") +  
  scale_x_discrete(labels=c("0" = "None", "1" = "One",
                              "2" = "Two")) +
  #scale_y_continuous(limits = quantile(d$fixDur, c(0.1, 0.9))) +
  stat_compare_means(comparisons = my_comparisons)+ # Add pairwise comparisons p-value
  stat_compare_means(label.y = 13) +     # Add global p-value
  geom_text(data = p_meds, aes(x = multNum, y = med, label = med), 
              size = 3, vjust = -1.5) +
  ggtitle("RT Distribution vs. Number of Multipliers: Absolute Summed Value < $0.50")

p0

RT by Mult for “Easy” trials

# create medians to insert as text
d<- S_M
d <- d[d$absSummedVal>1.0, ] # limit to absolute summed values below 0.50
x <- d
x$multNum = d$multNum+1
p_meds <- ddply(x, .(multNum), summarise, med = median(rt))
p_meds$med = round(p_meds$med, digits = 2)  # round to two decimal values

# List of conditions to compare
my_comparisons <- list( c("0", "1"), c("0", "2"), c("1", "2") )

p0 = ggboxplot(d, x = "multNum", y = "rt", color = "multNum", palette = "jco" )+ #outlier.shape=NA
  
  labs(y = "Total RT (seconds)", x = "Number of Multipliers") +
  theme(legend.position="none") +  
  scale_x_discrete(labels=c("0" = "None", "1" = "One",
                              "2" = "Two")) +
  #scale_y_continuous(limits = quantile(d$fixDur, c(0.1, 0.9))) +
  stat_compare_means(comparisons = my_comparisons)+ # Add pairwise comparisons p-value
  stat_compare_means(label.y = 13) +     # Add global p-value
  geom_text(data = p_meds, aes(x = multNum, y = med, label = med), 
              size = 3, vjust = -1.5) +
  ggtitle("RT Distribution vs. Number of Multipliers: Absolute Summed Value > $1.00")

p0

d<- S_M
d <- d[d$absSummedVal<0.25, ] # limit to absolute summed values below 0.50
subject_means <- group_by(d, multNum) %>%
  dplyr::summarize(median = mean(correct, na.rm = T))
subject_means

What is the mean abs value of combos with 0/1/2 mults?

Highest for 1 multiplier

d <- S_M

d0 <- d[d$multNum == 0, ]
d1 <- d[d$multNum == 1, ]
d2 <- d[d$multNum == 2, ]

mean(d0$absSummedVal)
mean(d1$absSummedVal)
mean(d2$absSummedVal)

PSYCHOMETRICS

Study 1/Study 2, RT and Accuracy

load("Data/NS_M.Rdata")

There were 50 or more warnings (use warnings() to see the first 50)

d1 <- NS_M
d2 <- S_M
# Create ID for each DF
d1$study <- "Standard Mult."
d2$study <- "Swap Mult."
# Need to uniquely number Subjects
d2$subject <- d2$subject + 100
# Concat DFs
common_cols <- intersect(colnames(d1), colnames(d2))
df = rbind(
  d1[, common_cols], 
  d2[, common_cols]
)
# GROUPBY
subject_means <- group_by(df, subject, study) %>%
  dplyr::summarize(accuracy = mean(correct, na.rm = T), rt = mean(rt, na.rm = T), meanSwap = mean(swapAmount, na.rm = T))

Error in summarise_impl(.data, dots) : 
  Evaluation error: object 'swapAmount' not found.

Create Binned Values to test against Accuracy and RT

setwd("/Users/djw/Dropbox/PHD/PRESENTATIONS/2017_SNE/Plots")

The working directory was changed to /Users/djw/Dropbox/PHD/PRESENTATIONS/2017_SNE/Plots inside a notebook chunk. The working directory will be reset when the chunk is finished running. Use the knitr root.dir option in the setup chunk to change the the working directory for notebook chunks.

ggsave("Prob.png", width = 19, height = 12, units = "cm")

Same as above but for abs Val for RT and Fixation

  
setwd("/Users/djw/Dropbox/PHD/PRESENTATIONS/2017_SNE/Plots")

The working directory was changed to /Users/djw/Dropbox/PHD/PRESENTATIONS/2017_SNE/Plots inside a notebook chunk. The working directory will be reset when the chunk is finished running. Use the knitr root.dir option in the setup chunk to change the the working directory for notebook chunks.

ggsave("RT.png", width = 19, height = 12, units = "cm")

FUNCTIONS for SUMMARY STATS

From: http://www.cookbook-r.com/Graphs/Plotting_means_and_error_bars_(ggplot2)/#Helper%20functions

## Norms the data within specified groups in a data frame; it normalizes each
## subject (identified by idvar) so that they have the same mean, within each group
## specified by betweenvars.
##   data: a data frame.
##   idvar: the name of a column that identifies each subject (or matched subjects)
##   measurevar: the name of a column that contains the variable to be summariezed
##   betweenvars: a vector containing names of columns that are between-subjects variables
##   na.rm: a boolean that indicates whether to ignore NA's
normDataWithin <- function(data=NULL, idvar, measurevar, betweenvars=NULL,
                           na.rm=FALSE, .drop=TRUE) {
    library(plyr)
    # Measure var on left, idvar + between vars on right of formula.
    data.subjMean <- ddply(data, c(idvar, betweenvars), .drop=.drop,
     .fun = function(xx, col, na.rm) {
        c(subjMean = mean(xx[,col], na.rm=na.rm))
      },
      measurevar,
      na.rm
    )
    # Put the subject means with original data
    data <- merge(data, data.subjMean)
    # Get the normalized data in a new column
    measureNormedVar <- paste(measurevar, "_norm", sep="")
    data[,measureNormedVar] <- data[,measurevar] - data[,"subjMean"] +
                               mean(data[,measurevar], na.rm=na.rm)
    # Remove this subject mean column
    data$subjMean <- NULL
    return(data)
}
## Summarizes data, handling within-subjects variables by removing inter-subject variability.
## It will still work if there are no within-S variables.
## Gives count, un-normed mean, normed mean (with same between-group mean),
##   standard deviation, standard error of the mean, and confidence interval.
## If there are within-subject variables, calculate adjusted values using method from Morey (2008).
##   data: a data frame.
##   measurevar: the name of a column that contains the variable to be summariezed
##   betweenvars: a vector containing names of columns that are between-subjects variables
##   withinvars: a vector containing names of columns that are within-subjects variables
##   idvar: the name of a column that identifies each subject (or matched subjects)
##   na.rm: a boolean that indicates whether to ignore NA's
##   conf.interval: the percent range of the confidence interval (default is 95%)
summarySEwithin <- function(data=NULL, measurevar, betweenvars=NULL, withinvars=NULL,
                            idvar=NULL, na.rm=FALSE, conf.interval=.95, .drop=TRUE) {
  # Ensure that the betweenvars and withinvars are factors
  factorvars <- vapply(data[, c(betweenvars, withinvars), drop=FALSE],
    FUN=is.factor, FUN.VALUE=logical(1))
  if (!all(factorvars)) {
    nonfactorvars <- names(factorvars)[!factorvars]
    message("Automatically converting the following non-factors to factors: ",
            paste(nonfactorvars, collapse = ", "))
    data[nonfactorvars] <- lapply(data[nonfactorvars], factor)
  }
  # Get the means from the un-normed data
  datac <- summarySE(data, measurevar, groupvars=c(betweenvars, withinvars),
                     na.rm=na.rm, conf.interval=conf.interval, .drop=.drop)
  # Drop all the unused columns (these will be calculated with normed data)
  datac$sd <- NULL
  datac$se <- NULL
  datac$ci <- NULL
  # Norm each subject's data
  ndata <- normDataWithin(data, idvar, measurevar, betweenvars, na.rm, .drop=.drop)
  # This is the name of the new column
  measurevar_n <- paste(measurevar, "_norm", sep="")
  # Collapse the normed data - now we can treat between and within vars the same
  ndatac <- summarySE(ndata, measurevar_n, groupvars=c(betweenvars, withinvars),
                      na.rm=na.rm, conf.interval=conf.interval, .drop=.drop)
  # Apply correction from Morey (2008) to the standard error and confidence interval
  #  Get the product of the number of conditions of within-S variables
  nWithinGroups    <- prod(vapply(ndatac[,withinvars, drop=FALSE], FUN=nlevels,
                           FUN.VALUE=numeric(1)))
  correctionFactor <- sqrt( nWithinGroups / (nWithinGroups-1) )
  # Apply the correction factor
  ndatac$sd <- ndatac$sd * correctionFactor
  ndatac$se <- ndatac$se * correctionFactor
  ndatac$ci <- ndatac$ci * correctionFactor
  # Combine the un-normed means with the normed results
  merge(datac, ndatac)
}
## Gives count, mean, standard deviation, standard error of the mean, and confidence interval (default 95%).
##   data: a data frame.
##   measurevar: the name of a column that contains the variable to be summariezed
##   groupvars: a vector containing names of columns that contain grouping variables
##   na.rm: a boolean that indicates whether to ignore NA's
##   conf.interval: the percent range of the confidence interval (default is 95%)
summarySE <- function(data=NULL, measurevar, groupvars=NULL, na.rm=FALSE,
                      conf.interval=.95, .drop=TRUE) {
    library(plyr)
    # New version of length which can handle NA's: if na.rm==T, don't count them
    length2 <- function (x, na.rm=FALSE) {
        if (na.rm) sum(!is.na(x))
        else       length(x)
    }
    # This does the summary. For each group's data frame, return a vector with
    # N, mean, and sd
    datac <- ddply(data, groupvars, .drop=.drop,
      .fun = function(xx, col) {
        c(N    = length2(xx[[col]], na.rm=na.rm),
          mean = mean   (xx[[col]], na.rm=na.rm),
          sd   = sd     (xx[[col]], na.rm=na.rm)
        )
      },
      measurevar
    )
    # Rename the "mean" column    
    datac <- rename(datac, c("mean" = measurevar))
    datac$se <- datac$sd / sqrt(datac$N)  # Calculate standard error of the mean
    # Confidence interval multiplier for standard error
    # Calculate t-statistic for confidence interval: 
    # e.g., if conf.interval is .95, use .975 (above/below), and use df=N-1
    ciMult <- qt(conf.interval/2 + .5, datac$N-1)
    datac$ci <- datac$se * ciMult
    return(datac)
}

Difficult, Very Difficult, Easy, Overall RT and Accuracy by MultNum

d <- S_M

There were 41 warnings (use warnings() to see them)

# Remove abs summed values >1.00 and <0.50
d <- d[(d$absSummedVal<=0.50) | (d$absSummedVal>=1.00),]
# Create Difficulty Level and Factor it
d$difficulty = 1  # easy level
d$difficulty[d$absSummedVal<0.5] = 2 # Difficult level
d$difficulty[d$absSummedVal<0.25] = 3 # V.Difficult level
# Factor conditions
d$multNum <- factor(d$multNum)
d$difficulty <- factor(d$difficulty)
# FOR T-TESTS
subject_means <- group_by(d, subject, difficulty, multNum) %>%
  dplyr::summarize(accuracy = mean(correct), rt = mean(rt))
subject_means

# Paired TTest
# Accuracy
# Easy vs Hard
mean(subject_means$accuracy[subject_means$difficulty == 1])

[1] 0.9519747

mean(subject_means$accuracy[subject_means$difficulty == 3])

[1] 0.6752014

sd(subject_means$accuracy[subject_means$difficulty == 1])

[1] 0.06967094

sd(subject_means$accuracy[subject_means$difficulty == 3])

[1] 0.1424739

       
t.test(subject_means$accuracy[subject_means$difficulty == 1],
       subject_means$accuracy[subject_means$difficulty == 3], paired = TRUE)

    Paired t-test

data:  subject_means$accuracy[subject_means$difficulty == 1] and subject_means$accuracy[subject_means$difficulty == 3]
t = 14.77, df = 68, p-value < 2.2e-16
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 0.2393794 0.3141673
sample estimates:
mean of the differences 
              0.2767733 

# Easy, 0 Mult/1 Mult
mean(subject_means$accuracy[subject_means$difficulty == 1 & subject_means$multNum == 0])

[1] 0.9877854

mean(subject_means$accuracy[subject_means$difficulty == 1 & subject_means$multNum == 1])

[1] 0.9482313

sd(subject_means$accuracy[subject_means$difficulty == 1 & subject_means$multNum == 0])

[1] 0.03119811

sd(subject_means$accuracy[subject_means$difficulty == 1 & subject_means$multNum == 1])

[1] 0.05554537

t.test(subject_means$accuracy[subject_means$difficulty == 1 & subject_means$multNum == 0],
       subject_means$accuracy[subject_means$difficulty == 1 & subject_means$multNum == 1], paired = TRUE)

    Paired t-test

data:  subject_means$accuracy[subject_means$difficulty == 1 & subject_means$multNum ==  and subject_means$accuracy[subject_means$difficulty == 1 & subject_means$multNum ==     0] and     1]
t = 5.1382, df = 22, p-value = 3.772e-05
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 0.02358934 0.05551889
sample estimates:
mean of the differences 
             0.03955411 

# Hard, 0 Mult/1 Mult
mean(subject_means$accuracy[subject_means$difficulty == 3 & subject_means$multNum == 0])

[1] 0.585889

mean(subject_means$accuracy[subject_means$difficulty == 3 & subject_means$multNum == 1])

[1] 0.7318523

sd(subject_means$accuracy[subject_means$difficulty == 3 & subject_means$multNum == 0])

[1] 0.1247233

sd(subject_means$accuracy[subject_means$difficulty == 3 & subject_means$multNum == 1])

[1] 0.1278067

t.test(subject_means$accuracy[subject_means$difficulty == 3 & subject_means$multNum == 0],
       subject_means$accuracy[subject_means$difficulty == 3 & subject_means$multNum == 1], paired = TRUE)

    Paired t-test

data:  subject_means$accuracy[subject_means$difficulty == 3 & subject_means$multNum ==  and subject_means$accuracy[subject_means$difficulty == 3 & subject_means$multNum ==     0] and     1]
t = -4.4311, df = 22, p-value = 0.0002105
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -0.21427733 -0.07764923
sample estimates:
mean of the differences 
             -0.1459633 

# Paired TTest
# RT
# Easy vs Hard
mean(subject_means$rt[subject_means$difficulty == 1])

[1] 2.714638

mean(subject_means$rt[subject_means$difficulty == 3])

[1] 3.428044

sd(subject_means$rt[subject_means$difficulty == 1])

[1] 0.9481159

sd(subject_means$rt[subject_means$difficulty == 3])

[1] 1.115059

       
t.test(subject_means$rt[subject_means$difficulty == 1],
       subject_means$rt[subject_means$difficulty == 3], paired = TRUE)

    Paired t-test

data:  subject_means$rt[subject_means$difficulty == 1] and subject_means$rt[subject_means$difficulty == 3]
t = -10.001, df = 68, p-value = 5.434e-15
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -0.8557461 -0.5710664
sample estimates:
mean of the differences 
             -0.7134063 

# Easy, 0 Mult/1 Mult
mean(subject_means$rt[subject_means$difficulty == 1 & subject_means$multNum == 0])

[1] 2.503676

mean(subject_means$rt[subject_means$difficulty == 1 & subject_means$multNum == 1])

[1] 2.666323

sd(subject_means$rt[subject_means$difficulty == 1 & subject_means$multNum == 0])

[1] 0.857097

sd(subject_means$rt[subject_means$difficulty == 1 & subject_means$multNum == 1])

[1] 0.8832649

t.test(subject_means$rt[subject_means$difficulty == 1 & subject_means$multNum == 0],
       subject_means$rt[subject_means$difficulty == 1 & subject_means$multNum == 1], paired = TRUE)

    Paired t-test

data:  subject_means$rt[subject_means$difficulty == 1 & subject_means$multNum ==  and subject_means$rt[subject_means$difficulty == 1 & subject_means$multNum ==     0] and     1]
t = -3.2763, df = 22, p-value = 0.003451
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -0.2656007 -0.0596928
sample estimates:
mean of the differences 
             -0.1626467 

# Hard, 0 Mult/1 Mult
mean(subject_means$rt[subject_means$difficulty == 3 & subject_means$multNum == 0])

[1] 3.485174

mean(subject_means$rt[subject_means$difficulty == 3 & subject_means$multNum == 1])

[1] 3.216965

sd(subject_means$rt[subject_means$difficulty == 3 & subject_means$multNum == 0])

[1] 1.173273

sd(subject_means$rt[subject_means$difficulty == 3 & subject_means$multNum == 1])

[1] 1.183658

t.test(subject_means$rt[subject_means$difficulty == 3 & subject_means$multNum == 0],
       subject_means$rt[subject_means$difficulty == 3 & subject_means$multNum == 1], paired = TRUE)

    Paired t-test

data:  subject_means$rt[subject_means$difficulty == 3 & subject_means$multNum ==  and subject_means$rt[subject_means$difficulty == 3 & subject_means$multNum ==     0] and     1]
t = 2.3206, df = 22, p-value = 0.02997
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 0.0285203 0.5078978
sample estimates:
mean of the differences 
              0.2682091 

#------------#
# PLOT       #
#------------#
# For RT
# Stats Summary
datac <- summarySEwithin(d, measurevar="rt", withinvars=c("multNum","difficulty"), idvar="subject")
ggplot(datac, aes(x=difficulty, y=rt, fill=multNum)) +
    geom_bar(position=position_dodge(.9), colour="black", stat="identity") +
    geom_errorbar(position=position_dodge(.9), width=.25, aes(ymin=rt-ci, ymax=rt+ci)) +
    coord_cartesian(ylim=c(2,4)) +
    labs(y = "Total RT (s)", x = "Difficulty (net value)") +
    scale_x_discrete(labels=c("1" = "Easy (>1)", "2" = "Difficult (0.5<0.25)",
                              "3" = "Very Difficult (<0.25)")) +
    scale_y_continuous(breaks=seq(2,4,0.2)) +
    theme_bw() +
    scale_fill_discrete(name="Number of\nMultipliers") 

    #ggtitle("RT vs. Difficulty + Number of Multipliers")  
# For Accuracy
# Stats Summary
datac <- summarySEwithin(d, measurevar="correct", withinvars=c("multNum","difficulty"), idvar="subject")
ggplot(datac, aes(x=difficulty, y=correct, fill=multNum)) +
    geom_bar(position=position_dodge(.9), colour="black", stat="identity") +
    geom_errorbar(position=position_dodge(.9), width=.25, aes(ymin=correct-ci, ymax=correct+ci)) +
    coord_cartesian(ylim=c(0.5,1)) +
    labs(y = "Accuracy", x = "Difficulty (net value)") +
    scale_x_discrete(labels=c("1" = "Easy (>1)", "2" = "Difficult (0.5<0.25)",
                              "3" = "Very Difficult (<0.25)")) +    
    scale_y_continuous(breaks=seq(0,1,0.1)) +
    theme_bw() +
    scale_fill_discrete(name="Number of\nMultipliers") 

    #ggtitle("Accuracy vs. Difficulty + Number of Multipliers")  

ANOVA

ANOVA on difference between multiplier and non multiplier trials RT

Based on this, there is no significant effect (as expressed through RT) in having one multiplier, however there is for having two.

Based on this, there is no significant effect (as expressed through RT) in having one multiplier, however there is for having two.

FIXATION DURATION

load("Data/S_M_K.Rdata")
d <- S_M_K

# Factor conditions
d$subject <- factor(d$subject)

# FOR T-TESTS
subject_means <- group_by(d, subject) %>%
  dplyr::summarize(firstFix = mean(fixDur[fixNum == 1]),
                   middleFix = mean(fixDur[fixNum > 1 & revFixNum > 1]),
                   finalFix = mean(fixDur[revFixNum == 1]))
subject_means

# Paired TTest
# RT
mean(subject_means$firstFix)
mean(subject_means$middleFix)
mean(subject_means$finalFix)
sd(subject_means$firstFix)
sd(subject_means$middleFix)
sd(subject_means$finalFix)

t.test(subject_means$middleFix,
       subject_means$finalFix, paired = TRUE)

Linear Models for Dependent (Fixed) Effects (not taking random effects into account)

First Fixation Duration vs. First Image Total Value

First Fixation Duration vs. First Mult

APRIL 24: NEW ANALYSES

Mixed Models

FOR STARTERS: Does summed value affect RT (using log RT)?

Perhaps unsurprisingly based on what we plotted before, summed value has a highly significant effect on reaction time (controlling for random effects of subjects).

Summed Value as the individual TOTAL values (the value x multiplier) of the FACE and HOUSE

Based on this there is a significant difference between mean RT and rt.model2 as well as between rt.model2 and rt.model1

What about interaction between faceTotal and houseTotal?

So there is a significant interaction bewteen the total house value and the total face value (as expected).

And then what if we look at the components (value * multiplier) of the Total Face and Total House Value?

Again, based on the anova analysis there seems to be significance in the interactions between the values and the multipliers.

—————————————–

RT PLOTS

Plot RT random effects of subjects

Plot RT with multiple lines (average of each subjects line)

• multiplier trials
• flip trials
• non flip trials
• non mult trials

Note that with the flip trials the summed values ranged from -0.6 to 0.6.

Plot Accuracy (% correct) with multiple lines (average of each subjects line)

• multiplier trials
• flip trials
• non flip trials
• non mult trials

Same Plot but just for “Difficult” choices

RT for difficult choices by multNum

CHOICE CURVE

% acceptance (sinusoid) -sinusoid for non mult/mult/flip

Value of first item vs. First item fixation time

• no mult
• 2x mult
• 3x mult

Value of second item vs. Second item fixation time

• no mult
• 2x mult
• 3x mult

Summed value vs second item fixation time

Number of swaps

• based on summed value

• based on ambiguity of individual stimuli (ie. closer to zero) -stimulus left has ambiguity X, stimulus right has ambiguity Y, summed value has ambiguity Z -how much does individual ambiguity vs combined ambiguity affect RT/swaps

SALIENCY TEST: Mean Fixation House v. Face

People looked at HOUSES longer. More salient? Or more ambiguous?

FINAL FIXATION: Tied to Value?

Are people taking longer for flip trials after accounting for fact that flip trials ALWAYS have multipliers (and non-flip trials don’t)

This is not currently working

Questionnaire Data

Plotting GPA vs. Earnings (unfiltered data):

Stroop Data

ANOVA for significance

#reformat Data Frame
found = which(rt_by_condition!=-999,arr.ind=T)
rtANOVA = data.frame(cbind(found,rt_by_condition[found]))
names(rtANOVA) = c('subj','cond','rt')

rtANOVA$subj = factor(rtANOVA$subj)
rtANOVA$cond = factor(rtANOVA$cond)

myaov = aov(rtANOVA$rt~rtANOVA$cond+Error(rtANOVA$subj))
summary(myaov)

T-Test

#reformat Data Frame
t.test(rt_by_condition$congruent,rt_by_condition$incongruent,paired=T,mu=0,alternative="two.sided",var.equal=T)

Based on Anova/TTest seems like there is a significant difference.

#Select dataframe to use
d <- Stroop.df.clean

#mean RT and Final earnings by subject
Stroop.performance <- group_by(d, subject) %>%
  dplyr::summarize(rt = mean(Response.rt, na.rm = T), accuracy = mean(Response.corr, na.rm = T))
Stroop.performance$performance = Stroop.performance$rt * 1/Stroop.performance$accuracy
#invert so bigger nunbers are better
Stroop.performance$performance = 1/Stroop.performance$performance

Stroop Performance vs. Task Performance

vs Earnings

vs Accuracy

vs RT

EFFECT SIZE

Cohen’s D

#Select dataframe to use
d <- total_M_clean3

#mean RT and Final earnings by subject
rt_mults <- group_by(d, subject) %>%
  dplyr::summarize(mult_0 = mean(RT[multNum==0]), mult_1 = mean(RT[multNum==1], na.rm = T), mult_2 = mean(RT[multNum==2], na.rm=T))
rt_mults

library(lsr)
cohensD(rt_mults$mult_0, rt_mults$mult_1)
cohensD(rt_mults$mult_0, rt_mults$mult_2)

MEANS

List of Mean RTs for Mults

#Select dataframe to use
d <- total_M_clean3

#mean RT depending on multiplier combination
for(i in 1:3){
  for(j in 1:3){
    m <- mean(d$RT[d$mult1House==i & d$mult2Face==j])
    cat(sprintf("House Mult = %s and Face Mult = %s\n", i, j))
    cat(sprintf("Mean: %f\n\n", m))
  }
}

MORE MIXED EFFECTS STUFF

Simulate RTs based on data

Abs Val vs. RT

Subject Level

APA Format Plotting

apatheme=theme_bw()+
  theme(panel.grid.major=element_blank(),
        panel.grid.minor=element_blank(),
        panel.border=element_blank(),
        axis.line=element_line(),
        text=element_text(family='Times'))

LME

with help from Liz

ggplot(fitDf, aes(`Fixed Effects`, z, fill=`Fixed Effects`)) + 
  geom_bar(stat = "identity", width = 0.5) + 
  geom_errorbar(aes(ymin=z-se, ymax=z+se), width=0.4) +
  geom_text(aes(label=star), colour="black", vjust=0, size=6) +
  scale_x_discrete(limits = positions) +
  theme_minimal() +
  theme(axis.title.x=element_text(size=14),
      axis.title.y = element_text(size = 14))+
  theme(legend.position="none")
setwd("/Users/djw/Dropbox/PHD/PRESENTATIONS/2017_SNE/Plots/")

The working directory was changed to /Users/djw/Dropbox/PHD/PRESENTATIONS/2017_SNE/Plots inside a notebook chunk. The working directory will be reset when the chunk is finished running. Use the knitr root.dir option in the setup chunk to change the the working directory for notebook chunks.

ggsave("fixedEffects.pdf", width = 20, height = 12, units = "cm")

Test to see if RT for incorrect is longer than RT for correct?

LOOKING AT QUESTIONNAIRE DATA

Note that this is the for the “cleaned” subjects. Running this with all of the subjects gives a significant effect to Self-Control.

MULT DIF and MODELS…

TEST TO SEE IF INTERACTION IS SIGNIFICANT

F-test is used to compare the residual sum of squares of both the models

drop1(model.1i, test = “F”) *doesn’t seem to work with lme (example is with lm)

MODEL VALIDATION

1. residuals versus fitted values to verify homogeneity

2. a QQ-plot or histogram of the residuals for normality

3. residuals versus each explanatory variable to check independence

**Instead of a visual inspection, it is also possible to apply a test for homogeneity. Sokal and Rohlf (1995) describe three such tests, namely 1. the Barlett’s test for homogeneity *sensitive to non-normality! 2. Hartley’s Fmax test and the log-anova 3. Scheffe ́-Box test

help(t.test)